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Abstract. We analyse the clustering features of Large Scale Structures (LSS) in the presence of 
massive neutrinos, employing a set of large-volume, high-resolution cosmological N-body simulations, 
where neutrinos are treated as a separate collisionless fluid. The volume of 8/i“^Gpc^, combined 
with a resolution of about 8 x 10^° h~^ Mq for the cold dark matter (GDM) component, represents 
a significant improvement over previous N-body simulations in massive neutrino cosmologies. In 
this work we focus, in the first place, on the analysis of nonlinear effects in GDM and neutrinos 
perturbations contributing to the total matter power spectrum. We show that most of the nonlinear 
evolution is generated exclusively by the GDM component. We therefore compare mildly nonlinear 
predictions from Eulerian Perturbation Theory (PT), and fully nonlinear prescriptions (halofit) 
with the measurements obtained from the simulations. We find that accounting only for the nonlinear 
evolution of the GDM power spectrum allows to recover the total matter power spectrum with the 
same accuracy as the massless case. Indeed, we show that, the most recent version of the halofit 
formula calibrated on AGDM simulations can be applied directly to the linear GDM power spectrum 
without requiring additional fitting parameters in the massive case. As a second step, we study the 
abundance and clustering properties of GDM halos, confirming that, in massive neutrino cosmologies, 
the proper definition of the halo bias should be made with respect to the cold rather than the total 
matter distribution, as recently shown in the literature. Here we extend these results to the redshift 
space, finding that, when accounting for massive neutrinos, an improper definition of the linear bias 
can lead to a systematic error of about 1-2% in the determination of the linear growth rate from 
anisotropic clustering. This result is quite important if we consider that future spectroscopic galaxy 
surveys, as e.g. Euclid, are expected to measure the linear growth-rate with statistical errors less than 
about 3% at z < 1. 
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1 Introduction 

In the ACDM cosmological model the three active neutrinos of the standard model of particle physics 
are assumed to be massless. Nevertheless, already in 1998 the Super-Kamiokande collaboration pre¬ 
sented evidence of neutrino oscillations [1], indicating that at least two neutrinos are massive, and, 
more recently, new neutrino oscillation experiments seem to exclude a vanishing flavor mixing angle at 
more than 10 cr (see, e.g. [2, 3]). Indeed, the study of the effects of neutrino masses on cosmological 
observables is of particular relevance for two, distinct, reasons. First, the absolute neutrino mass 
scale remains unknown and, in this respect, cosmology plays a key role in its determination, being 
gravity sensitive to the total neutrino mass, Sm^, rather than to the mass splitting. Second, an 
accurate description of massive neutrino effects on LSS is required to avoid systematic errors in the 
determination of cosmological parameters, as the dark energy density and equation of state, whose 
measurements represent one of the main goals of current and future cosmological experiments. 

Massive neutrino cosmologies have been extensively studied in the literature (see [4, 5] for a 
review). In particular the linear perturbation theory in the presence of massive neutrinos is well 
understood and it is widely used to derive constraints on E and other cosmological parameters, 
from present and future CMB and galaxy surveys observations [6-23]. However, the increasing pre¬ 
cision of cosmological parameter measurements further requires an accurate description of nonlinear 
corrections, which can be obtained by means of a direct analysis of the output of N-body simulations 
accounting for a massive neutrino component. This kind of simulations turns to be computationally 
very expensive and quite challenging if a good mass resolution and a large box size are required at the 
same time in order to build realistic mock catalogues for present and future galaxy and weak lensing 
surveys. 

The “Dark Energy and Massive Neutrino Universe” (DEMNUni) simulation project presented in 
this and in the companion paper [24] , addresses in a consistent way the problem of structure formation 
in massive neutrino cosmologies, and represents the state of the art of neutrino simulations in terms 
of volume and mass resolution. At present, as explained in more details in §3 below, these simulations 
include a baseline ACDM model and three cosmologies characterised by different values of Em,y, 
with all the cosmologies sharing the same total matter density 12^, as well as the same amplitude 
of primordial curvature perturbations. In the near future, the DEMNUni set will be extended with 
the inclusion of an evolving dark energy background, with different equations of state rc, in order to 
study the degeneracy between the neutrino mass S rrii, and w at the nonlinear level. 
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This work focuses on the analysis of the clustering features of LSS in massive neutrino cos¬ 
mologies, as extracted from the DEMNUni simulation set, including, in the first place, an accurate 
comparison of the nonlinear total matter power spectrum with the most recent predictions in PT, over 
the range of scales where such predictions are expected to significantly improve over linear theory. To 
this end we will separately consider the clustering of the CDM and neutrino components as well as 
their cross-correlation, following previous indications that such distinctions are qnite relevant for an 
accurate description of gravitational instability in these cosmological models (see, e.g. [25-27]). At the 
same time, we will test different implementations of popular fitting formulae such as halofit [28-30] 
describing the nonlinear regime. A specific study of non-Gaussian aspects of the matter distribution 
in massive neutrinos cosmologies is left to future work [31]. 

In addition, we will revise previous results on the abundance and clustering of dark matter halos, 
taking advantage of the large volume and mass resolution provided by the DEMNUni simulations. 
In this respect, our results provide further evidence that a description in terms of cold (rather than 
total) dark matter perturbations is indeed required for the extension of common fitting functions 
for the halo mass function to the case of massive neutrinos. In this respect, despite the significant 
systematic effect that the wrong assumption would lead to the interpretation of cluster abundance 
observations [32] , these early indications have been so far ignored in the analysis of relevant data-sets 
(e.g. [33]). We believe that, whether or not particle-based neutrino simulations provide an faithful 
description of reality, such systematic uncertainties should be at the very least accounted for in any 
exhaustive error budget. Similar considerations can be extended to the interpretation of galaxy power 
spectrum measurements where the assumption of a constant linear bias with respect to the CDM or 
total matter power spectrum can lead to quite different constraints on the total neutrino mass [26]. 
In this work, for the first time we extend from real to redshift space our results on halo clustering in 
the presence of massive neutrinos, providing a first insight into possible, additional systematic effects 
due to erroneous assumptions in the modelling of galaxy clustering anisotropies in massive neutrino 
scenarios. 

This paper is organised as follow. In §2 we give a brief description of linear perturbation theory in 
the presence of massive neutrinos. In §3 we present the N-body runs and the details of the numerical 
setup. In §4 we present our measurements of all the components to the matter power spectrum, along 
with comparison with PT results and popular fitting formulae. In §5 and §6, we focus our attention 
on the mass function and the power spectrum of halos both in real and redshift space. Finally, in 
§7 we summarise the main results of the present work and draw our conclusions, highlighting several 
applications of the DEMNUni simulations to the testing of different observables and probes in the 
context of massive neutrinos cosmologies. 

2 Massive neutrinos perturbations 

Massive neutrinos are considered hot dark matter, since they decouple in the early Universe as rela¬ 
tivistic particles, just before the onset of Big Bang Nucleosynthesis, and become non-relativistic only 
at later times. This transition happens when the neutrino temperature drops below the value of their 
mass at redshift [34] 

1 + Znr^ 1890 . (2.1) 

After this time they contribute to the total energy density of the Universe as dark matter so that the 
total matter energy density is given by 

Ilm — Hcrfm 4“ Hfe 4“ Hi/ : (2.2) 

where, in addition to the CDM and baryon components, accounts for the energy density associated 
to the massive neutrino component. Since we will consider scales much larger than the Jeans length 
of baryons, in what follows we will define Sic = ^cdm 4- Slti as the density corresponding to the sum 
of CDM and baryon densities, and, for simplicity, will generically refer to it as a single “CDM” 
contribution, denoted by the subscript “c”. It can be shown that is proportional to the total 
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neutrino mass as 


(2.3) 


= 




93.14/i2eV ’ 

where the proportionality factor depends on the assumed photon temperature and neutrino to photon 
temperature ratio and it should be evaluated numerically in the most general case [35]. As a result, 
at late times, i.e. 2 <C z„j., the effect of neutrinos on the expansion rate of the Universe is completely 
degenerate with a change of the CDM and baryon components. 

On the other hand, at the perturbation level massive neutrinos have a peculiar effect on matter 
density fluctuations. Defining the density contrasts for neutrino and CDM respectively as <5^ = Spujpn 
and 5c = 5pc/pc, where the total mass density is pm = Pm + 5pc + 5pv, pm = pc + pv being the total 
background matter density, we can write 


(5m = (1 - fu)5c + fu5u 


(2.4) 


where fu represents the neutrino fraction defined as It follows that the total matter 

power spectrum can be written as the sum of three contributions corresponding respectively to the 
CDM power spectrum, Pec, the neutrinos power spectrum, and the cross-power spectrum between 
CDM and neutrinos, Pcu{5c5l), that is 


= (1 - /,)2 Pec + 2U (1 - U) Pc. + fl P. 


(2.5) 


which shows that the neutrino fraction has a direct impact on the total matter power spectrum by 
modifying the relative contributions of the two components. Eqs. (2.2) to (2.5) introduce the notation 
adopted throughout the rest of the paper. 

Over the age of the Universe, neutrinos travel an average distance that depends on their thermal 
velocity and, in turn, on their mass. This free streaming length determines the scale below which 
neutrinos density perturbations are washed-out, and is given by (see, e.g. [5]) 


Afs(wi.,z) ~ 8.1 


Ho{lPz) (leY 
H{z) 


1,-1 


Mpe. 


( 2 . 6 ) 


Notice that, for particles becoming non-relativistic during matter domination, as it is usually 
the case for neutrinos, the comoving free streaming length, Apg/a, is actually decreasing in time, 
and therefore assumes the largest value at the time of the non-relativistic transition. This peculiar 
distance corresponds to the wave-number 


/1 V \ 

knr = kFs{Znr) ^ 0.018 j /lMpc“^ . 

\m. J 


(2.7) 


This scale is typically larger than the scale at which nonlinear effects manifest themselves at low 
redshifts. At any redshift, scales larger than 1/knr are affected by the presence of massive neutrinos, 
and, in first approximation, we can write 


Pmmik^ — 


\ Peeik) 


for k kr. 


^ (1 - /i.)^ Pcc{k) for fc > , 


( 2 . 8 ) 


while the exact value of the damping scale will retain a residual redshift dependence. In the left 
panel of Figure 1 we plot the ratio Pmm/Pcc for Y,mv= 0.3, 0.53 eV at redshifts z = 0, 2, showing 
the two asymptotic regimes of Eq. (2.8). On very large scales the ratio goes to one, while at small 
scales it approaches (1 — regardless of the redshift. Intermediate scales are instead affected by 
the actual value of the free streaming wave number kpsim., z) = 27ra/AFs(w;y, z), and by its redshift 
dependence. 

It can be shown [36] that, with respect to the massless case, the total linear power spectrum, 
Pmm, in massive neutrino scenarios is asymptotically suppressed at z = 0 by a constant factor on 
scales k knr 


Pmm {k', f.) 
Pmm{k-, f. = 0) 


- 1 - 8/, 


(2.9) 
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Figure 1. Linear theory results in massive nentrino cosmologies. Left panel. Ratio of the total matter power 
spectrnm to the CDM power spectrum at redshifts z = 0 {continuous curves) and 2 = 2 {dashed curves) for 
two different valnes of the sum of neutrino masses, E 0.3 eV in red and E 0.53 eV in green. Dotted 
lines denote the asymptotic value at small scales of (1 — fu)^. Right panel ratio at 2 = 0 of the total matter 
power spectrum {continuous curves) and CDM power spectrum {dashed curves) for the same two cosmologies 
to the ACDM prediction. 


while from Eq. (2.8) and Eq. (2.5), it follows that the suppression for the CDM power spectrum, 
is given by a factor ^ (1 — Qf^)- The difference in the suppression between the two power spectra is 
shown in the right panel of Figure 1. 

3 Simulations 

The DEMNUni simulations have been conceived for the testing of different probes, including galaxy 
surveys, CMB lensing, and their cross-correlations, in the presence of massive neutrinos. To this 
aim, this set of simulations is characterised by a volume big enough to include the very large-scale 
perturbation modes, and, at the same time, by a good mass resolution to investigate small-scales 
nonlinearity and neutrino free streaming. Moreover, for the accurate reconstruction of the light-cone 
back to the starting redshift of the simulations, it has been used an output-time spacing small enough 
that possible systematic errors, due to the interpolation between neighbouring redshifts along the line 
of sight, result to be negligible. 

The simulations have been performed using the tree particle mesh-smoothed particle hydrody¬ 
namics (TreePM-SPH) code GADGET-3, an improved version of the code described in [37], specifically 
modified in [38] to account for the presence of massive neutrinos. This version of GADGET-3 follows 
the evolution of CDM and neutrino particles, treating them as two distinct collisionless fluids. For the 
specific case of the DEMNUni simulations, a GADGET-3 version, modified for OpenMP parallelism 
and for memory efficiency, has been used to smoothly run on the BG/Q Fermi cluster. 

Given the relatively high velocity dispersion, neutrinos have a characteristic clustering scale 
larger than the GDM one. This allows to save computational time by neglecting the calculation of the 
short-range tree-force induced by the neutrino component. This results in a different scale resolution 
for the two components, which for neutrinos is fixed by the PM grid (chosen with a number of cells 
eight times larger than the number of particles), while for CDM particles is larger and given by the 
tree-force (for more details see [38] ). This choice does not affect the scales we are interested in, since, 
as shown in [39], the application of the short-range tree force to neutrino particles is required only 
to describe the neutrino density profile inside massive halos at low redshifts. In this work we are not 
considering this kind of analysis on such small scales. 
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Y,m^[eV] 

^cdm 


^8,mm 

*^S,cc 

m^p[h ^ Mq] 

mp[h ^Mq] 

SI 

0.00 

0.2700 

0.000 

0.846 

0.846 

8.27 X 10^° 

— 

S2 

0.17 

0.2659 

0.013 

0.803 

0.813 

8.16 X 10^° 

1.05 X 10^ 

S3 

0.30 

0.2628 

0.022 

0.770 

0.786 

8.08 X 10^° 

1.85 X 10^ 

S4 

0.53 

0.2573 

0.040 

0.717 

0.740 

7.94 X 10^° 

2.28 X lO'^ 


Table 1. Summary of cosmological parameters and derived quantities for the four models assumed for the 
DEMNUni simulations. The values fit = 0.05, Qm = 0.32, h = 0.67, Us = 0.96 are shared by all the models. 


Each DEMNUni run starts at redshift Zin = 99, and is characterised by a comoving volume of 
8 h~^ Gpc^, filled with 2048^ dark matter particles and, when present, 2048^ neutrino particles. Given 
the large amount of memory required by the simulations, baryon physics is not included. Nevertheless, 
since we are looking for deviations from a fiducial reference ACDM model, it is expected that baryon 
feedback cancels out when considering relative effects, i.e. the massive with respect to the massless 
neutrino case. In addition we expect also that any additional effect produced by the interplay of 
neutrinos with baryon physics should be of higher order. This is supported also from [40] which shows 
that the neutrino induced suppression in the total matter power spectrum is very much the same also 
when neutrinos are considered in the presence of baryons. Therefore, our choice should not affect the 
results presented in this work. 

We have produced a total of four different DEMNUni simulations, choosing the cosmological 
parameters according to the Planck 2013 results [41], namely a flat ACDM model generalised to a 
z/ACDM , le. a massive neutrino model, by varying only the sum of the neutrino masses over the 
values = 0, 0.17, 0.3, 0.53 eV (and consequently the corresponding values of Dy and Dc, while 
keeping fixed Dm and the amplitude of primordial curvature perturbations Ag). Table 1 provides 
a summary of the cosmological parameters that characterise the different runs. It shows as well 
the derived quantities as,mm and as,cc corresponding to the r.m.s. of perturbations on spheres of 
8h“^Mpc, computed respectively for the total and CDM matter components, and the CDM and 
neutrino mass particle resolutions, which vary according to the value of and 

The simulations are characterised by a softening length e = 20 h~^ kpc, and have been run on the 
Eermi IBM BG/Q supercomputer at CINECA'^, employing about 1 Million cpu-hrs per simulation 
(including the production of halo and sub-halo catalogues). For each simulation, 62 outputs have 
been produced, logarithmically equispaced in the scale factor a = 1/(1 -I- 2 ), in the redshift interval 
z = 0 — 99, 49 of which lay between z = 0 and z = 10. For each of the 62 output times, it has been 
produced on-the-fly one particle snapshot, composed by both CDM and neutrino particles, one 3D 
Cartesian grid of the gravitational potential, (p, and one 3D cartesian grid of its time derivative, (p, 
with a mesh of dimension 4096^ that covers a comoving volume of 2 h~^ Gpc, for a total of about 
90 TB of data per simulation. In addition, in order to build parent-halo catalogues, the simulation 
outputs have been first processed with the friends-of-friends (FoF) algorithm included in the GADGET- 
3 package. The code is applied only to CDM particles with linking length set to 0.2 times the mean 
inter-particle distance. A minimum number of 32 particles to identify a structure has been assumed, 
fixing the minimum parent-halo mass to MpoP — 2.5 x 10^'^ h~^ Mq. The FoF catalogues have 
subsequently been processed with the subfind algorithm [42, 43], which identifies locally overdense, 
gravitationally bound regions within an input parent halo. With this procedure some of the initial 
FoF parent halos are split into multiple sub-halos, with the result of an increase of the total number 
of identified objects, and of a lower minimum mass limit. Here, a minimum number of 20 particles 
has been adopted in order to constitute a valid sub-halo. Moreover, by means of a specific routine 
included in the subfind algorithm, spherical overdensity halo catalogues have been produced, and, in 
particular, for each halo the value of R 200 has been computed; this is defined as the radius enclosing 
a matter (CDM-|-neutrinos) density equal to 200 times the mean density of the Universe Pm{z) at 

^Consorzio Interuniversitario del Nord-Est per il Calcolo Automatico, Bologna, Italy: 
http: // www.hpc. cineca. it / services / iscra. 
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Figure 2. Comparison between the DEMNUni runs and previous, recent simulations of massive neutrino 
cosmologies in terms of CDM particle mass resolution and simulation volume. Grey diagonal lines indicate 
the number of CDM particles. 


redshift z, and the corresponding mass in terms of M 200 = 200 pm ^ ^ 200 • We stress here that all 
the post-processing has been modified to take the neutrino component properly into account. 

Figure 2 presents a comparison of the DEMNUni simulations to previous N-body simulations 
with massive neutrino particles, in terms of the CDM particle mass resolution nip and simulation 
volume V. Diagonal grey lines indicate the number of CDM particles. With respect to previous 
simulations, the DEMNUni suite represents an improvement of about an order of magnitude in terms 
of particle number (only [44] considered a larger box, but with considerably smaller mass resolution). 

4 Matter power spectra 

As shown in §2, the shape of the linear power spectrum is quite sensitive to the value of neutrino 
masses, and this dependence becomes even stronger in the mildly and fully nonlinear regimes [40, 44]. 
Taking advantage of the large DEMNUni simulations volume, in this section we aim at testing the 
accuracy of current analytical predictions for the nonlinear matter power spectrum, Pmnn in the 
presence of massive neutrinos. To this end, we measure individually the different components to Pmm 
in Eq. (2.5), from very large scales, k ^ 0.003 h Mpc“^, down to fully nonlinear scales, k Mpc“^, 
and compare these measurements with PT predictions, in the mildly nonlinear regime, and fitting 
functions as halofit [30, 45], in the fully nonlinear regime^. The goal here is to understand if possible 
departures from the linear regime of neutrino perturbations have to be taken into account for precision 
cosmology at the % level. 

Before proceeding to the comparison with the nonlinear, analytical predictions, however, we take 
a look at each component and its relative contribution to the nonlinear Pmm measured from the 
simulations. Figures 3 and 4 show the CDM auto power spectrum, (1 — P^c {dashed curves), the 
neutrino auto power spectrum, {dotted curves), and the cross CDM-neutrino power spectrum, 

2 /j, (1 — /j,) Pm, {dot-dashed), as extracted from the simulations {thick curves), and the corresponding 

^While the mass resolution of the DEMNUni simulations would allow to look at much smaller scales, of the order of 
k ~ 10/iMpc“^, we do not investigate this regime since it is dominated by baryon physics [46—50]. 
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k [h Mpc ’] k [h Mpc *] 


Figure 3. Comparison between the different contributions to the nonlinear matter power spectrum, (1 — 
Pcc {dashed curves), Pw (dotted) and 2 f,^ (1 — /i,) P^u (dot-dashed), as described in the text. All the 
measurements at redshifts 2 = 0, 1, 2 are shown with shades varying, respectively, from blue to red. Thin 
coloured curves correspond to the respective linear predictions. Dashed and dotted grey lines on the top 
panels show the shot-noise contributions to the CDM and neutrinos power spectra. The shaded area in the 
bottom panel shows values below the 1-cr, relative, Gaussian uncertainty on AP(k)/P(k) = l/y^27rfc^/fcp kf 
being the fundamental frequency of the simulation box. 


linear predictions (thin curves). Each plot corresponds to a different value of the neutrino mass in the 
simulations. Within each plot, the top panel shows the adimensional power spectrum, dir P(k), 
the middle panel the ratio of each contribution to Pmm, in the nonlinear and linear cases, and the 
bottom panel shows the ratio (Pnl — Pl)/ Pmrm be. the difference between each measured nonlinear 
component, Pml^ in Eq. (2.5), and the corresponding linear prediction, P^, compared to the nonlinear 
Pmm- All the measurements at redshifts z = 0, 1, 2 are shown with different shades varying from blue 
to red. Dashed and dotted grey lines on the top panels show the shot-noise contributions to the CDM 
and neutrinos power spectra, respectively. The shaded area in the bottom panel show the region 
below the 1 -ct, Gaussian uncertainty on P(k), given by A^P(fc) = P^(fc)/(27rfc^/fc^), kf being the 
fundamental frequency of the simulation box (the shot-noise contribution to the variance is ignored 
for simplicity). 

In the linear regime, fc < 0.1 ZiMpc”^, most of the power comes from the CDM component, with 
the cross power spectrum term accounting roughly for a 10% of the total Pmm-, and the neutrinos 
being already negligible (see middle panels of Figures 3 and 4). On nonlinear scales, from the lower 
panels we notice that most of the nonlinear contribution to Pmm is still given by Pcc, which is the 
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k [h Mpc ’] k [h Mpc *] 


Figure 4. Same as figure 3, but for 0.3, 0.53 eV. 


only component significantly deviating from the linear prediction, while the other two terms remain 
linear, and, therefore, less and less important as we move to smaller scales^. It should be kept in mind 
that each contribution in Figures 3 and 4 is multiplied by the proper power of fi, or (1 — f^), with 
/y ~ 0{2 — 4%) . We therefore conclude that, on the scales probed by the DEMNUni simulations, 
and considering the present constraints on the sum of neutrino masses, the total nonlinear matter 
power spectrum, in massive neutrino cosmologies, can be described at the 1 % level by accounting for 
the nonlinear evolution of CDM perturbations alone, while adopting the linear prediction for the other 
components. This result will be useful for analytical predictions of the nonlinear Pmm discussed in 
the next sections. 

4.1 Perturbation Theory 

For the first time in the literature, the large volume of the DEMNUni simulations gives us the 
possibility to measure the matter power spectrum in neutrino cosmologies at the 1% accuracy level, 
on a very large range of scales, allowing a test of PT predictions at the accuracy level required by 
current and futures galaxy surveys. 

®Let us notice that measurements of the neutrino auto power spectrum and neutrino-cold matter cross power 
spectrum shown in both Figures 3 and 4, present a spurious contribution at small scales not to be confused with a 
residual shot noise component. Such contribution, which scales roughly (but not exactly) as 1/fc, is due to the neutrino 
velocity distribution and it does not grow with time, so that at low redshift is overtaken by the expected neutrino power 
spectrum. A detailed discussion of numerical issues associated with the set-up of the initial conditions will be discussed 
in a forthcoming work [51]. 
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Several works in the literature have discussed the effects of massive neutrinos in cosmological 
perturbation theory beyond the linear level [10, 27, 52-58]. In these descriptions, the neutrino compo¬ 
nent is treated, similarly to the CDM one, as a single perfect fluid, fully characterised in terms of its 
density and velocity divergence (see, however, [59, 60] for a multiple-flow approach to the evolution 
of neutrino perturbations). The main difference with respect to the ACDM case is represented by 
an effective sound speed modifying the Euler equation for the neutrino component, and accounting 
for the neutrino velocity distribution. The first consequence, at the linear level, is a scale-dependent 
linear growth factor, D{k, z), for both the CDM and neutrino components. However, the perfect-fluid 
approximation fails to provide an accuracy for the neutrino power spectrum below the 10% level [57]. 
Nonetheless, as shown in §4, since the neutrino contribution to the total matter power spectrum is 
order of magnitudes smaller than the CDM one, such discrepancies on the neutrino component alone 
do not affect significantly the CDM and total matter power spectra. Therefore, we will assume the 
two-fluid approximation for all the comparisons of analytical versus numerical results in this section. 

In addition, even if in the mildly nonlinear regime the effective sound speed affects as well the 
mode-coupling at all the orders of the perturbative expansion, we will follow the same approach 
adopted by [54]. They have shown that limiting the neutrino-induced scale-dependence to the linear 
growth factor alone (and, therefore, the use of standard EdS- like kernels in the perturbative expan¬ 
sion) proves to be a quite good approximation to the full PT solution for the nonlinear CDM field, 
on scales where PT is expected to be accurate. 

Finally, we will make the additional approximation, already proposed in [52], of describing the 
neutrino perturbations by means of their linear solution. While this is not per se a good assumption 
[27, 58], it does provide the correct neutrino contribution to the total matter power spectrum on the 
(large) scales where such contribution is relevant. 

As a starting point for future, more accurate comparisons, we will consider, therefore, the fol¬ 
lowing perturbative prediction for the total matter power spectrum 

P^lik) = (1 - Ur pPJik) + 2 (1 - U) U PtAk) + fl PAik) ■ (4.1) 

Here, the contribution PAA{k,z) represents the nonlinear power spectrum predicted in perturbation 
theory along the lines of [54], i.e. it is computed in terms of its linear counterpart, PA{k,z), which 
provides the correct linear scale-dependence of the growth factor, but assumes the standard EdS 
nonlinear kernels in the perturbative expansion. Differently from previous works, however, we do 
not only consider standard, one-loop corrections to P^, but we take into account also standard PT 
two-loop corrections, as well as the “regularised” predictions based on the multi-point propagator 
expansion of [61], computed using the RegPT code of [62]. Since the RegPT code does not account 
for the evolution of the scale-dependent linear growth in massive neutrino models, we produce all the 
predictions at z > 0 by providing the corresponding z-input linear power spectrum as a “fake z = 0 
input” required by the code. 

While this is not the most rigorous approach, it represents a practical application, to massive 
neutrino scenarios, of available tools developed within the ACDM framework. As we will see, the 
gain in accuracy achieved by recent resummation schemes, applied here to the CDM component 
alone, might compensate for the crude approximations that this approach implies. Clearly we are 
only considering predictions for the CDM and total matter power spectra, as these statistics are the 
relevant ones for galaxy clustering and weak leasing observations. 

In Figure 5 we show the perturbative results against the measurements at z = 0.5, 1, and for 
T,m^= 0, 0.53 eV. Error bars are the theoretical expectation for a Gaussian field, that is 


AP^{k) 


1 

2Tr k'^ kf 


P{k) + 



(4.2) 


where kf = 21: jL is the fundamental frequency of the simulation box, L being its linear size, and n 
is the particle number density accounting for the shot-noise component^. 

^The relatively small scatter of data points with respect to the error bars is due to the specific seed chosen for the 
random number generator used for the set-up of the initial conditions [63]. 
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Figure 5. Perturbation Theory predictions for the cold matter power spectrum Pcc{k). Each panel shows the 
measurements from the N-body simulations, divided by the reference power spectrum given by the two-loop, 
standard PT results {black points with error-bars). Also shown are the corresponding ratios for the linear 
{green, dotted), one-loop, standard PT {blue, thin, dashed curve), multi-point propagator expansion at one- 
{red, thick, dashed) and two-loops {red, thick, continuous) ass obtained from the RegPT code of [62]. 


Let us first notice that, in the ACDM case (left panels in Figure 5), the two-loops standard 
PT does not provide a good fit to the data at low redshifts [64], while it reproduces fairly well the 
simulation measurements at z > 1. Analytic predictions are 1% accurate at z = 1, up to a maximum 
wave-number k^ax — 0.3/iMpc~^ 

Turning our attention to the CDM power spectrum in massive neutrino cosmologies (right panels 
in Figure 5), at z > 1 we find approximately that the 1% accuracy is recovered up to a value of k^ax 
very close to the ACDM case. On the other hand, at lower redshifts, z < 1, we notice that, with 
respect to the massless case, PT predictions are in better agreement with the measurements extracted 
from the DEMNUni simulations. This is a consequence of the fact that, due to massive neutrino free 
streaming, the nonlinear evolution of CDM perturbations is suppressed, and therefore the range of 
scales in which PT corrections show the same accuracy as in standard ACDM cosmologies increases 
with f,y, i.e. with the relative contribution of the total neutrino mass to Qm- In addition, since, massive 
neutrino free streaming has an impact already at the linear level, these models are characterised by a 
lower value of Ug (and erg than a ACDM universe with the same 

5 It is worth noticing that the agreement may also depend on the simulation mass resolution; we expect that much 
higher resolutions lead to more power at small scales [65]. 
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Figure 6. Same as Figure 5 but for the total matter power spectrum with Em^= 0.3 eV {left panels) and 
0.53 eV {right panels). Theoretical predictions have been obtained using Eq. (4.1), i.e. computing nonlinear 
correction for CDM only. 


Given a prediction for the CDM power spectrum accurate at the 1% level up to a given kmax, we 
check if the perturpative approach of Eq. (4.1) reproduces, with the same level of accuracy, the total 
matter power spectrum measured from the simulations. The results are illustrated in Figure 6, which 
shows that, indeed, the linear treatment of the Pd, and P^, contributions to the total Pmm proves to be 
a very good approximation. The difference in the accuracy of the predictions between the E m^= 0.3 
eV {left panels) and the E?7i^= 0.53 eV {right panels) cases is again mainly due to the different values 
of /i, and, therefore, to the different effect of neutrino free streaming on dark matter perturbations, 
according to the total neutrino mass: for a given value of a larger value of E not only increases 
the relative amount of neutrino perturbations that are washed out below the free streaming scale, 
ApS) consequently reducing the contribution of P^, and P^, to the total Pmrm but also decreases the 
factor (1 — /j^)^ in front of Pec in Eq. (4.1), where the nonlinear evolution of Pec is in turn suppressed, 
with respect to the massless case, by the action of the total gravitational potential sourced both by 
CDM and massive neutrinos. Apart from inducing a scale-dependence of the linear growth factors for 
CDM and total matter, the main direct product of this physical mechanism is represented by a lower 
amplitude of linear perturbations at z = 0, where Us.cc = 0.786, 0.740, and ers.mm = 0.770, 0.717, for 
E m,j= 0.3, 0.53 eV, respectively. We will show in §4.2 that, on scales much smaller than the so-called 
turn-over scale, beyond the mildly nonlinear regime, k > 0.2 /iMpc”^, where the growth factor scale- 
dependence induced by neutrino free streaming approaches its asymptotic value depending only on 
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/y (see Eq. (6.12) in §6.2), the effect of massive neutrinos on Pec and Pmm mostly reduces to a mere 
rescaling of the power spectra in the massless case, according to the values of us.cc and CTg.mm- 

Recently the BOSS collaboration released new constraints on neutrino masses based on measure¬ 
ments of the galaxy power spectrum multipoles in redshift space [21] at the mean redshift z = 0.57. 
The BOSS analysis was based on a prediction for the total matter power spectrum computed applying 
standard PT directly to the linear Pmm- The outcome of such a calculation, in principle even less 
theoretically justified than our crude assumptions, is very similar to our results obtained via Eq. (4.1); 
we have checked that any difference between the two approaches stays below the 1% level at the scales 
relevant for current observations. 

A test somehow similar to the one presented here is shown in [66], where the authors compare 
different PT predictions, including the Time-RG method of [55, 67], to simulations of CDM particles, 
modifying only the background evolution and the initial conditions to account for free streaming 
massive neutrinos. They show that PT predictions are in agreement with the measurements of CDM 
power spectra, extracted from their CMD simulations, at the % level, over a similar range of scales 
as tested in this work. However they assume a scale-independent growth rate to rescale back the late 
time (z = 0) CDM power spectrum, Pcc(fc), to the initial redshift of the simulations. By doing so, the 
linear dynamics cannot be correctly recovered at any z other than z = 0. 

Finally, the crucial results of this section rely on the discussion of §4 and the measurements 
shown in Figures 3 and 4, that is the contributions from the two terms, other than Pec, entering 
Eq. (4.1), and in particular from the cross power spectrum P^y, always remain negligible compared 
to Pec on nonlinear scales, at least for the level of accuracy requested for PT to be useful. 

4.2 Fitting functions 

Given the limitations of the perturbative approaches, it is sometimes convenient and/or sufficient to 
resort to fitting functions for the nonlinear matter power spectrum. In this section we see how the 
approximation of linear evolution for neutrino perturbation can be applied as well to the halofit 
prescription [45]. 

The HALOFIT formula, originally based on stable clustering considerations [68, 69], provides a 
mapping between the linear power spectrum and the nonlinear one, which in turn depends on few 
cosmological parameters, e.g., and several free parameters determined by comparisons against 
measurements from N-body simulations. A new, more accurate, version of the fitting formula has 
been recently presented by [30]. The revised formula is expected to be accurate at the 5% level for 
k < l/iMpc”^ and z < 10, while it degrades to the 10% level for k < lO/iMpc”^ and z < 3. 

In the context of massive neutrinos cosmologies. Ref. [40] provided corrections (controlled by a 
few, additional parameters) to the original halofit formula, in order to account for neutrinos effects 
on the nonlinear total matter power spectrum®. However, as shown in Figures 3 and 4, nonlinear 
corrections to the cross power spectrum between CDM and neutrinos are below the percent level, 
therefore we wonder if, similarly to Eq. (4.1), a fitting formula for the total matter power spectrum, 
Pmm{k), could be given in terms of the linear predictions, Pcy{k) and P^j^{k), and the halofit fitting 
formula applied directly to the linear CDM power spectrum, P^^{k), that is 

P^^{k) = (1 - UY Pg^{k) + 2 /. (1 - U) Pi{k) + fl P^k) . (4.3) 

Here P^^(fc) = pHF[Pcc{k)], where the halofit mapping, Phf, between linear and nonlinear power 
spectra does not account for additional corrections due to massive neutrino physics. 

In the left column of Figure 7 we show the ratio of the measured CDM power spectrum, Pcc{k), 
to the prediction, P^^(k). Each row corresponds to a distinct value of Em,/, while each panel shows 
the value of this ratio at redshifts z = 0, 0.5, 1, 1.5, 2, with colour shades ranging from blue to red. 
The shaded areas denote the regions beyond the accuracy claimed for the revised formula of [30]. 
The left panels of Figure 7 show, indeed, not only that P^^{k) provides the expected accuracy, but 
that it works equally well for all the considered values of the total neutrino mass. This is essential 

®The massive neutrino corrections of [40] have been implemented, along the revised version of [30], in the most recent 
versions of the camb [70] and class [71] codes. 
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Figure 7. Left column-, ratio of the measured CDM nonlinear power spectrum, Pec, to the halofit prediction 
Pcc,HF from [30], with no additional parameters to account for neutrino effects. Each panel correspond to one 
value of T,mu, with different colours denoting different redshifts: 2 = 0, 0.5, 1, 1.5, 2, blue to red. Shaded 
areas denote the regions beyond the accuracy expected for the formula. Right column-, ratio of the measured 
total matter nonlinear power spectrum, Pmm, to the prediction, P^mi of Eq. (4.3). 


to justify our assumption of applying the halofit mapping to the CDM component alone. Here, we 
stress again that the version of halofit employed for J^hf[PccW] does not include any effect due to 
massive neutrinos on the CDM clustering, since here we use the halofit version developed by [30] 
in the ACDM framework. This result is similar to that obtained in §4.1 for perturbation theory, and 
it is a crucial step before checking the validity of the assumptions made in Eq. (4.3). 

The right column of Figure 7 presents the ratio of the measured total matter power spectrum, 
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Figure 8. Comparison between the nonlinear power spectra measured in two ACDM models {blue curves) 
with (78 = 0.68 (continuous) and erg = 0.70 (dashed) against the same quantity measured in massive neutrinos 
models (red curves) with ag^mm = 0.68 (continuous) and (Jg,cc = 0.70 (dashed). All the quantities are shown 
as ratios to the power spectrum for a ACDM model with ag = 0.83 to compare the shape of the suppression 
resulting from either a lower overall normalisation or neutrino effects. Measurements for this figure only are 
from the simulations of [72]. 


Pmm{k), to the prediction P^^(k) of Eq. (4.3). In addition, dashed curves show the inverse ratio 
of P^^(k) to the specific halofit extension to massive neutrino cosmologies of [40], denoted as 
Pmm~‘'ik)- We notice that the simple prescription of Eq. (4.3), while avoiding introducing additional 
parameters to the fitting formula of [30] , remains well within the expected halofit accuracy. On the 
other hand, the discrepancies between the prediction of Eq. (4.3) and the P^^~^(k) fit of [40] are 
within 4%, with the latter ad-hoc fit not improving particularly over the former. 

Here we would like to make some final considerations. The results of this section rely on the fact 
that, in the first place, in the absence of strong baryon physics, on scales much beyond the mildly 
nonlinear regime, k > 0.2 hMpe”^, where the linear growth factor scale-dependence induced by 
neutrino free streaming approaches its asymptotic value, the extent of the nonlinear evolution of the 
power spectrum can be accounted for via the amplitude of the linear power spectrum, regardless of the 
physical mechanism responsible for the amplitude itself. In other terms, on nonlinear scales, we expect 
a similar behaviour for the nonlinear matter power spectrum evolved from a linear power spectrum 
suppressed by massive neutrino free streaming, as from a different linear power spectrum with a lower, 
primordial normalisation. This assertion can be easily tested with N-body simulations and the results 
are shown in Figure 8, which makes use, however, of measurements from the simulations described 
in [72, 73]. The plot shows ratios, with respect to the matter power spectrum in a given ACDM 
cosmology with erg = 0.83, of the matter power spectrum in other ACDM cosmologies (with different 
(Tg normalisations ), and of the CDM and total matter power spectra in massive neutrino cosmologies, 
where either the value of erg^cc or ug^mm are matched to the ACDM ones. We observe that on nonlinear 
scales, 0.2 < k < 1 h Mpc~^, a ACDM model with a given erg is nearly indistinguishable from a massive 
neutrino model with the same value for erg cc- A lower agreement is found when the match is done 
in terms of as,mm, *-e. in terms of the normalisation of the linear total matter power spectrum, since 
the relevant nonlinear evolution is given by CDM perturbations, and the amplitude normalisation can 
not take correctly into account the contribution from remaining linear neutrino perturbations (this 
should be taken into consideration especially in weak-lensing analyses, which probe directly the total 
matter power spectrum). This result represents the well known degeneracy between Y,mi, and Ug 
at small scales, and implies that the enhanced, nonlinear suppression of the matter power spectrum 
on nonlinear scales can be obtained tuning the normalisation of the linear one, without resorting to 
peculiar effects of massive neutrinos. Nonetheless, here we stress that, such kind of degeneracy can 
be broken when observing the power spectra on a much larger range of scales, 0.01 < k < 5 /iMpc”^, 
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by means of future large sky galaxy surveys as, e.g., Euclid. 


5 Halos mass function 


The abundance of massive clusters as a function of redshift is a key cosmological probe for neutrino 
mass measurements [74-76]. In particular, massive neutrinos could alleviate the tension between the 
value of (Tg inferred from the primary CMB anisotropies [41, 77] and the one obtained from low redshift 
probes [18, 19, 32, 33, 78-81]. A careful calibration of the halo mass function in simulations with a 
massive neutrino component is, therefore, required to avoid systematic effects in the determination of 
cosmological parameters. 

The overall larger volume and higher mass resolution of the DEMNUni simulations with respect 
to previous studies (see Figure 2) allow an improved analysis of the abundance of massive clusters 
over several decades in mass. In this work, we identify halos in two different ways: using a Friends-of- 
Friends (FoF) algorithm with linking length equal to 0.2 times the mean inter-particle distance, and 
with the SUBFIND-SO code included in the gadget package, which identifies halos whose spherical 
overdensity is 200 times the mean background matter density, M 200 = (47r/3) i ?200 200 p, as already 
mentioned in §3. FoF halo masses have been corrected following the empirical prescription by [82], 
and, for the present analysis, no halos with less than 50 particles have been considered. 

The halo mass function, i.e. the number density n{M) of halos of mass between M and M + dM, 
is often expressed as 


n{M) 


Pr. 

M 


/(o-. z) 


diner ^ 
dM 


dM , 


(5.1) 


where most of the cosmological information is encoded in the variance of the matter distribution in 
the linear regime 

a^{R,z)= j d^kP{k,z)Wl{k), (5.2) 


smoothed on the scale R with a Top-Hat filter in real-space W/j(fc)^. The function f{a,z) can be 
either predicted in the Press-Schechter framework [83] (see [84, 85] for recent results) or fitted to 
numerical simulations (see, e.g. [86-88]). The smoothing scale and the halo mass M are related by 
the choice of the filter W{kR), and, in the Top-Hat case, is given by 



(5.3) 


Let us notice that in Eq. (5.1) as in Eq. (5.3), in the case of massive neutrino models, we need to use 
the background density of the cold (rather than total) matter component pc to define halo masses, 
since, as largely shown in the literature, the contribution of neutrinos (both bounded and unbounded) 
to the mass of CDM halos is completely negligible [25, 44, 89]. 

Another potential ambiguity in massive neutrino scenarios is that the variance cr(i?, z) in Eq. (5.1) 
could in principle be computed for cold or total matter perturbations. However, as shown in [25, 26], 
the halo mass function is well described by analytic predictions, and by common fitting functions 
developed in the context of ACDM cosmologies, only if cr = CTcc, *-e. the variance of the matter 
distribution is computed from the cold matter linear power spectrum Pcc{k). Here we summarise the 
physical mechanism behind this choice, referring the reader to [25, 26, 90, 91] for further details. In 
the spherical collapse model, assuming general relativity, the evolution of a region is controlled by the 
amount of mass inside its initial volume. For a Gaussian random field, as the linear density field d, a 
natural choice to describe the system is therefore its variance smoothed on the scale R associated to 
the initial region, a{R). Due to the tiny value of neutrino masses, the scale R is much smaller than the 
free streaming length of massive neutrinos Apsj typically tens of Megaparsec, and therefore neutrino 
perturbations are vanishingly small inside the collapsing region. This implies that cr(i?) should be 
computed using the CDM plus baryon perturbations only. A similar argument applies to the critical 

^To partially remove volume effects we use the fundamental frequency of the box, kf = 27r/2000hMpc“^, as the 
lower integration limit in Eq. (5.2). 
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Figure 9. Massive neutrino effects on the mass function of FoF halos. Data points show the mass function 
measured from the DEMNUni simulations for = 0, 0.17, 0.3, 0.53 eV. Errors are derived from the 

assumption of a Poisson distribution for each bin. Theoretical predictions are obtained in terms of the fitting 
formula of [87] as a function of acc{M) {continuous curves) and amm{M) {dashed curves). Different panels 
correspond to redshifts ^ = 0.5, 1. 


overdensity required for collapse at a given redshift, Scr, that depends on neutrino masses only through 
their effect on the background evolution [25]. Note that CTcc > o-mm always (see Figure 1) implying 
that for a given cosmology, using a^c in the computation of the mass function leads to predicting more 
halos than using amm (see [92] for implications in galaxy cluster counts observations). Of particular 
importance is that, as shown in [26], the universality of the function /(cr, z) with respect to different 
cosmologies is recovered only for cr = acc- 

In Figure 9 we show the suppression of the FoF halo abundance due to free streaming massive 
neutrinos. The points show measurements from the DEMNUni simulations, at two redshifts z = 0.5, 1, 
while the theoretical predictions are obtained in terms of the fitting formula of [87] as a function of 
acciM) {continuous curves) and amm{M) {dashed curves). Let us notice that the volume and mass 
resolution of the DEMNUni simulations allow to recover the mass function of the FoF halos on a quite 
large range of scales, 5 x 10^^ ^ M < 10^^ M©. 

As a further inspection, in Figure 10 we show the ratios, over the ACDM case, of the FoF halo 
mass function measured from the DEMNUni simulations at four different redshifts, z = 0, 0.5, 1, 1.5. 
Again, we compare these measurements to the predictions from the MICE fit to FoE halos in [87]. 
Note that the MICE fit has been obtained using simulations of ACDM cosmologies alone, and it is 
non-universal in redshift, i.e. f{o',z) does present an explicit redshift dependence on top of the cr(z) 
one. 

For the massive neutrino cosmologies we find that the cold matter prescription, with p = and 
O’ = CTcc in Eq. (5.1), reproduces very well the measurements when the MICE formula is assumed 
for the function /(tr, z). In fact, we find the same level of agreement encountered for the ACDM 
model, therefore showing the same level of non-universality with respect to redshift. At high masses 
the measurements in the simulations drop below the MICE fit, because of the smaller volume of 
the DEMNUni simulations with respect to the MICE simulations. Our mass function measurements 
largely confirm the findings of [26], but we stress that the DEMNUni simulations allow to extend such 
analysis by a factor of 8 in resolution while accounting, due to their larger volume, for the large-scale 
modes particularly relevant for this kind of studies [87]. 

We repeat the analysis for Spherical Overdensity (SO) halos, identifying halos as spherical regions 
with a mean matter density equal to two hundred times the background density, as described in §3. 
The measurements, shown as ratios to the ACDM case, are presented in Eigure 11, where now we use 
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Figure 10. Massive neutrino effects on the mass function of FoF halos. Data points show the ratio of the 
mass function measured from the DEMNUni suite for Emiy= 0.17, 0.3, 0.53 to the massless case. Errors are 
derived from the assumption of a Poisson distribution for each bin. Theoretical predictions are obtained in 
terms of the fitting formula of [87] as a function of a^ciM) {continuous curves) and {dashed curves). 

Different panels correspond to redshifts z = 0, 0.5, 1, 1.5 {clockwise from the top-right panel). 


the Tinker et al. fitting formula as a reference for SO halos [86]. Also in this case we observe that 
the Tinker et al. fit, developed in the ACDM framework, provides a good fit to simulations with a 
massive neutrino component when p = and a = cr^c are used instead of their total dark matter 
counterparts. Finally, we can assert that a parametrisation of the halo mass function in terms of the 
CDM power spectrum gives a much better description of measurements from simulations, including 
massive neutrinos, with respect to the DM power spectrum, both for FoF and SO halos, therefore 
confirming the results of [25, 26, 93]. 

6 Halo bias and redshift space distortions 

Measurements of the linear bias of FoF halos in massive neutrinos simulations have been presented 
already in [26], where the authors showed how halo bias at large-scales is scale-independent, as in 
standard ACDM cosmologies, only if defined with respect to the cold, rather than total, matter 
perturbations. It follows that the halo bias defined with respect to the total matter fluctuations 
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Figure 11. Same as Figure 9 but for Spherical Overdensity halos. Theoretical predictions, using a^c 
{continuous curves) and using a mm {dashed curves), have been obtained using the formula of [86]. 


presents a spurious scale-dependence simply arising from the difference between the cold and total 
matter power spectra, see Fig. 1. 

In this section we check, in the first place, these previous results taking advantage of the higher 
mass resolution of the DEMNUni simulations, with the aim of analysing halo population of sensibly 
lower mass. In addition, we perform a preliminary analysis of Redshift Space Distortions (RSD) 
in massive neutrino models with a specific attention to possible consequences of the aforementioned 
ambiguity in the halo bias definition. More specifically, we test the Kaiser limit [94] of RSD, as 
the simplest, observable feature of anisotropic galaxy clustering directly dependent on the proper 
definition of bias. The main goal is to highlight possible systematic effects in the determination of 
the growth rate, which, in the presence of massive neutrinos, is not simply a function of time, as 
in ACDM models, but presents a mild scale-dependence as well. This analysis is therefore rather 
complementary to the one of [89] where systematic errors on the determination of RSD parameters 
induced by neglecting neutrino masses have also been studied. 

While we have analysed, for this purpose, three distinct halo populations characterised by distinct 
mass thresholds, here we focus on the lower one, corresponding to M > 10^^ h~^ Mq, as less affected 
by shot-noise and characterised by smaller nonlinear bias corrections. We consider, however, this 
threshold at different redshifts, presenting the measurements at z = 0.5 and z = 1 as the most 
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Figure 12. Linear halo bias from the halo auto-power spectrum as a function of scale for halos of mass M > 
10^^ h~^ Mq in the DEMNUni simulations at redshift 2 ; = 0.5 {left panel) and 2 = 1 {right panel). Continuous 
curves with error bars show the bias defined w.r.t. the cold matter power spectrum as b{k) = Phh{k)/Pcc{k), 
while dashed curves denote the bias defined w.r.t. the total matter perturbations, b{k) = s/Phh{k)/Pmm{k)■ 
Horizontal lines show the value of the linear bias at fc ~ 0.06/i Mpc“^, assumed here, approximately, as the 
limit of validity for linear theory. 


relevant for future observations. 

6.1 The halo power spectrum in real space 

In the local bias approach [95], one can expand the halo overdensity Sh as a power series in Sm, such 
that on very large scales the bias relation can be written as 


Sh^bdjn, ( 6 . 1 ) 

where the constant b is the linear bias factor for the halo population. From the equation above (and 
from the definition of power spectrum) it follows, in ACDM cosmologies, a linear relation between 
the matter and halo power spectra given by 

Phh{k)-b'^Pmm{k). (6.2) 

Linear bias can therefore be measured in N-body simulations as b^ = Phh/Pmm in terms of the 
measured power spectra Phh and Pmm- In massive neutrino cosmologies, however, the P{k) at the 
RHS of the above expression can be either the cold or total matter power spectrum, leading to two 
possible definitions of linear halo (and galaxy) bias, that is 

‘-^1) (0,3) 

or 

= (6.4) 

V ^mm 

Analogous definitions can be given in terms of the cross halo-matter power spectra, bc^x = Phc/Pcc 
and bm,x = Phm/Pmm- The results from §5 and an analogy with ACDM cosmologies suggest that 
bias coefficients are scale-independent only if defined with respect to the cold matter field. 

In Figure 12 we present the measurements of halo bias on linear scales from the DEMNUni 
simulations for halos of masses M > h~^ Mq at redshifts 2 = 0.5, 1. We consider the bias 
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determined according to both Eq. (6.3) and Eq. (6.4). In all the cases, the bias defined with respect 
to the total matter perturbations exhibits a larger scale-dependence than the one defined w.r.t. the 
cold matter component, confirming the findings in [26]. In the rest of the section we will explore 
possible implications of this fact for the description of RSD at large scales. 


6.2 The halo power spectrum in redshift-space 

The biasing between the galaxy and matter distributions is not the only effect to be taken into 
account for a correct estimate of the matter power spectrum from galaxy redshift observations. In a 
real survey, the proper motions of galaxies with respect to the homogeneous expansion of the Universe 
affect the determination of their distance along the line-of-sight (see, e.g. [94, 96-98]). However, on 
sufficiently large scales these redshift space distortions provide information on the peculiar velocity 
field of matter perturbations, in particular on the growth rate of the density field, and, as a result, 
they are extensively used to constrain cosmological parameters and test deviation from standard 
gravity [99-102]. In preparation for future large spectroscopic surveys, the effect of neutrino masses 
on RSD modelling needs to be carefully investigated to avoid fake signatures of Modified Gravity. In 
the following we present a preliminary assessment of the scale-dependence of the growth rate induced 
by the free streaming length of neutrinos, and its measurement from the DEMNUni simulations. 

In the standard cosmological model, and in the large-scale limit where linear theory applies, the 
distortion induced by peculiar velocities on the density contrast Sm can be written in Fourier Space 
as [94] 

Sm,sik) = {1 + ffi^)Smik) , (6.5) 

where ^ = k ■ z/k is the angle between the line-of-sight and the wave vector k, while /(a) is the linear 
growth rate, defined as the logarithmic derivative of the linear growth factor D{a), that is 


/(«) 


d\nD{a) 
din a 


( 6 . 6 ) 


A similar expression holds for linearly biased tracers which follow the matter flow (be. with no velocity 
bias®). For instance, the halo overdensity in redshift space, Sh,s, can be written as 


Sh,s{k) = ib + ffi'^)Sm = {1 + I3fi^)bdm{k), 


(6.7) 


where b is the linear scale-independent bias, and we define /? = f /b. It follows that the halo power 
spectrum in redshift space can be written as the product of a polynomial in fi times the halo power 
spectrum in real space. Since the halo power spectrum in real space depends only on the modulus 
of the wave vector, the angular dependence induced by RSD is entirely encoded in the pre-factor of 
Eq. (6.7). As a result, the redshift-space halo power spectrum can be decomposed in multipoles using 
just the first three even Legendre polynomials Lf(/r) 


PhhAk) = {^ + Pk^''?Phh{k)= Y. Phh,iLiA)- ( 6 . 8 ) 

/^ 0 , 2,4 


The monopole, quadrupole and hexadecapole coefficients read 


PhhAk) = {l + \p+\A^Phh{k) 

PhhAk) = Phh{k) 


PhhAik) = —APhhik), 
35 


(6.9) 

( 6 . 10 ) 

( 6 . 11 ) 


^Note that the assumption of no velocity bias is better justified for halos w.r.t. the cold matter perturbations rather 
than w.r.t. the total matter ones. 
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Figure 13. Comparison of the ratio between the monopole of the redshift-space halo power spectrum to 
the real-space halo power spectrum, Phh,o{k)/Phh{k) {continuous curves) against the value predicted by the 
Kaiser formula, Eq. (6.9), as a function of /3 = /c/&c {dashed curves) or /3 = fmlbm {dotted curves). fc{k) and 
fm{k) correspond to the linear theory growth rate respectively of cold and total matter perturbations, while 
the be and bm are the measured values of the halo bias according to the two definitions of Eqs. (6.3) and (6.4). 
Notice that we do not use best fit values for the bias parameters, but we use instead the measured h{k) as a 
function of scale. 


where the RSD parameter /3 determines the relative amplitude of the multipoles. In a ACDM cos¬ 
mology, /3 is predicted to be scale-independent on large scales, as a direct consequence of the scale- 
independence of both linear bias and growth rate. On the other hand, both quantities might be 
scale-dependent when massive neutrinos are present. 

In the first place, in massive neutrino scenarios, the growth rate of cold matter perturbations, 
/c, becomes scale-dependent, as a natural consequence of the scale-dependent growth function D^. In 
particular, the small-scale asymptotic suppression expected in linear theory is given by [103] 

1 J'g _ ^25-24/,) - 1 - ^ , (6.12) 

and corresponds to a 2.4% effect for the largest value of the neutrino fraction assumed here, that is 
~ 0.039 in the 0.53 eV model. However, for the same model, the suppression is below the 

percent level on large scales, k < 0.03/iMpc“^. In the case of the total matter growth rate, fm{k), 
the suppression, again for the Sm^= 0.53 eV model, reaches the 1% level at slightly smaller scales, 
k ~ 0.05 Mpc“^ (linear theory predictions for both the ratio / c//acdm and fm/J acdm are shown 
as dashed curves in Figure 14). 

In the second place, if we defined the halo bias w.r.t. the total matter field, according to Eq. (6.4), 
we would add to the quantity f3 an additional source of scale-dependence. Let us notice here that, while 
one would consider reasonable that the choice of be, Eq. (6.3), is consisted with the choice of / = /c 
in the definition of /3 = f /b, we cannot exclude a priori other possibilities (including, for instance, a 
mixed (3 = fm/bc), since the physical origins of the two contributions in Eq. (6.7) are in fact distinct. 
However, in what follows, we will limit ourselves to compare the two cases {3 = fc/bc and (3 = fm/bm, 
the latter corresponding to naively neglecting any distinctions between matter contributions. 

In all the cases, we expect the Kaiser factor /3 to exhibit a scale-dependence, which, by choosing 
/3 = fc/bc, will result only from the growth rate fc{k). In the other case, (3 = fm/bm, a partial 
compensation occurs between the linear bias, which increases with the wave-number, and the growth 
factor that actually decreases with k. In Eigure 13 we show the ratio of the monopole of the redshift- 
space halo power spectrum, Phhp(k), to the halo power spectrum in real space, Phhik), for the 
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Figure 14. Top panels: Linear growth rate obtained from measnrements of the redshift-space (monopole) 
and real-space halo power spectrum according to Eq. (6.9) and from measurements of the linear halo bias 
b{k) = \/Phh{k) / Pcc{k) at redshifts 2 = 0.5 {left) and 2 = 1 ( right). Dashed cnrves denote the linear theory 
(scale-dependent) predictions for cold matter perturbations, fc{k). All the quantities are shown in terms of 
their ratio to the corresponding ACDM valnes. Bottom panels: Same as in the top panels, but assuming the 
linear halo bias defined in terms of total matter perturbations. In this case, the comparison is made to the 
linear theory prediction for the total matter growth rate, fm{k). 


models listed in Table 1, and we compare the measurements from the DEMNUni simulations to the 
corresponding predictions provided by the Kaiser formula written both in terms of /3 = fc/bc and 
/? = fm/bm- In such predictions we do not use best fit values for the bias parameters, but we use 
instead the measured b{k) where the scale-dependence is affected by cosmic variance. While the noisy 
measurements do not allow to clearly determine the scale-independence of such a quantity, we notice 
that /3 = fc/bc does provide a slightly better prediction than {3 = fm/bm, as compared to simulations. 

As a further test, in Figure 14 we consider directly the growth rate / obtained from measurements 
of the monopole of the redshift-space halo power spectrum, Phh,o{k), and the real-space halo and cold 
matter power spectra, Phh{k) and Pcc{k) respectively, under the assumption that the Kaiser formula 
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provides a good description of RSD at large scales. Assuming /3 = //5c, from Eq. (6.9) we have 


f{k) = 


jPhhik) 1 

/dS20 5 

/ Peeik) 3 

[V Phhik) \ 


(6.13) 


while assuming /3 = f/b^ we obtain a similar expression where Pmm replaces Pec- In order to 
further reduce cosmic variance, we show the recovered value of / as ratio to the ACDM case. These 
measurements, as a function of fc, can then be directly compared with linear theory predictions for 
fc{k) and /m(fc). 

The top panels of Figure 14 show / obtained from Eq. (6.13) under the cold matter hypothesis 
compared to the predictions for /c(fc) at redshift z = 0.5, 1. We notice that, in the higher redshift 
case in particular, despite the noise, the measurements are consistent with the linear predictions at 
large scales when / is obtained from the measured Pec- -A. greater discrepancy instead is observed 
when / is obtained from the measured total matter Pmm(^), as compared to the predictions for fm{k) 
in the bottom panels. If such discrepancy will be confirmed by future investigations, we could expect 
that using bm as a definition for bias in the Kaiser formula could lead to a systematic error on the 
determination of the growth rate at the level of 1-2%. 

Clearly this represents a very preliminary analysis of possible systematic effects in the determi¬ 
nation of the growth rate f(k) in the context of massive neutrinos cosmologies. If the description 
of halo clustering in terms of cold matter perturbation does indeed represents the correct approach, 
this test can serve as a confirmation that the DEMNUni simulations reproduce linear theory predic¬ 
tions for the growth rate including, to a certain extent, its scale-dependence. We reserve for future 
work a more detailed analysis of RSD effects on the matter density field as well as on realistic mock 
galaxy distributions, in massive neutrino scenarios. Of particular interest would be any description 
of nonlinearity that could extend theoretical predictions to smaller scales. 

7 Conclusions 

In this work we have presented an analysis of the LSS clustering features extracted from the DEMNUni 
simulation suite, a set of large-volume and high-resolution N-body simulations that include massive 
neutrinos as a particle component. The volume and resolution of these simulations allow the test of 
several cosmological observables, from the galaxy power spectrum to weak leasing statistics and CMB 
secondary anisotropies [24]. In this respect, the DEMNUni simulations represent the largest effort to 
date to include the effects of neutrinos mass in numerical predictions of the large-scale structure. As 
such, they serve, in the first place, the purpose of investigating our ability to constrain the properties 
of such important particle for both the standard cosmological model as the standard model of particle 
physics. At the same time, they can help assessing the effects of massive neutrino on the expected 
accuracy and precision in the determination of other cosmological parameters, first and foremost those 
describing possible departures from the ACDM paradigm. 

We have analysed basic quantities such as the matter power spectrum, the halo mass function, 
and halo clustering both in real and redshift space, comparing them with the simplest and most 
common theoretical predictions. This analysis is intended as a first step toward a more accurate 
test of the latter, taking advantage, however, of recent developments in the understanding of halo 
formation and clustering. 

The most fundamental quantity, in this respect, is the matter power spectrum. In Section 4 
we presented measurements of the power spectrum of each individual component of the total matter 
density, that is the cold dark matter and baryons (here treated on the same footing and called CDM) 
and the neutrinos, along with their cross-power spectrum and the total matter power spectrum, given 
by the weighted sum of the three. The analysis in terms of distinct components allows to test early 
hypothesis regarding the possibility of neglecting the nonlinear evolution of neutrino perturbations in 
analytical predictions [52]. In this respect we point-out, with direct reference to our measurements, 
that the neutrino power spectrum, as well as the cross-power spectrum provide a significant contribu¬ 
tion to the total matter power spectrum, Eq. (2.5), only at large-scales, where the linear approximation 
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is sufficient. Since the total matter power spectrum, Pmmik), and the cold matter power spectrum, 
Pcc{k), are the only quantities directly related to actual observables, a 1% accuracy appears to be 
achievable even neglecting the nonlinear evolution of neutrino perturbations. With this in mind, we 
focused on the accurate description of nonlinearity of the cold matter component alone, considering 
various predictions in perturbation theory, and using, as input quantity, the simple cold matter linear 
power spectrum. While this is clearly an effective, not rigorous, approach, we have shown that the 
accuracy provided by PT techniques, and their limits of validity on the matter power spectrum in 
ACDM cosmologies, can be achieved as well for massive neutrino models in a rather simple way. In a 
completely similar fashion, fitting formulae obtained from ACDM simulations, and used beyond the 
perturbative regime like halofit, can be safely applied to predictions of the nonlinear cold matter 
power spectrum, and simply extended, with the addition of the linear contributions from the neu¬ 
trino auto- and CDM-neutrino cross-spectra, to the total matter power spectrum, retaining the same 
accuracy in presence of massive neutrinos as in a standard massless cosmology, without resorting to 
additional fitting parameters. 

The analysis of the DEMNUni simulations has allowed us to confirm the results of [26, 93] on 
halo abundances in massive neutrino cosmologies, i.e. that the halo mass function is better described 
in terms of the variance of cold matter fluctuations alone. Similarly, from the study of halo clustering, 
we have confirmed as well that a constant linear, halo bias is recovered in the large-scale limit only 
when the bias relation is considered between halo and cold (rather than total) matter distributions. 
Such peculiar ambiguity in the definition of halo (and galaxy) bias, motivated a preliminary test of the 
Kaiser formula for redshift-space distortions [94], where linear bias is a relevant parameter. We found, 
despite the yet large statistical uncertainty, some indications of the expected scale-dependence of the 
growth rate of matter perturbation, characteristic of massive neutrino cosmologies. If such indications 
are correct, then we have shown that defining the galaxy bias w.r.t. the cold matter component is 
necessary to avoid a systematic error on the determination of the growth rate f{k) at the few percent 
level. We will return on the subject in more detail in future work. 
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